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In the 1950's Wigner and collaborators provided an explanation for the spectrum of hadronic excitations in 
terms of Random Matrix Theory. In the 1980's it was understood that some hadronic spectral properties were 
generic to systems whose classical limit is chaotic. We use string theory to demonstrate explicitly how, under 
very general conditions, recent holographic models of strong interactions have a spectrum compatible with 
Wigner's conjecture. 



I. INTRODUCTION 



During the 1950s Eugene Wigner set out to describe the general properties of energy levels of highly excited states of heavy 
nuclei [ 1 ]. The main idea was to describe such a complex nuclear system by a Hermitian Hamiltonian H and to connect the 
results to Random Matrix Theory (RMT) Q. Wigner proposed regarding a specific H as behaving like a large random matrix 
that is a member of an ensemble of Hamiltonians, all possessing similar general properties. Consequently, the spacings between 
energy levels of heavy nuclei could be more easily modelled by the spacings between successive eigenvalues of a random 
N x -matrix as N — » oo. 

During the 1980's a connection between RMT evolving from Wigner's work and quantum chaos was established through 
the understanding of simple and universal properties of the energy level fluctuations. Spectral fluctuations of quantum systems 
whose classical limit is fully chaotic show a strong level repulsion that depends only on the symmetries of the system, and 
agrees with the predictions of RMT Q. In contrast, classically integrable systems give rise to uncorrelated adjacent energy 
levels, that are well described by Poisson statistics (4). Haq, Pandey and Bohigas (HPB) demonstrated that the fluctuation 
properties of nuclear energy levels are in agreement with RMT [5 ]. HPB considered data consisting of 1407 resonance energies 
corresponding to 30 sequences of 27 different nuclei |5 ]; the connection with quantum chaos was then explored |6|. A statistical 
analysis looked at the experimentally measured mass spectrum of hadrons (TV, A, A and E and all the mesons up to /2(2340)) 
taken from the Particle Data Group Summary Tables and concluded that the nearest-neighbor mass-spacing distribution of the 
meson and baryon spectrum is described by the Wigner surmise corresponding to the Gaussian Othogonal Ensemble (GOE) (7). 
Lattice studies of QCD and its supersymmetric versions found similar results based on eigenvalues of the Dirac operator 0[9]]. 

Here, by studying the spectrum of hadronic states in holographic models developed in string theory, we provide an explicit 
realization of Wigner's conjecture. 

We first discuss the essential ingredients of the gauge/gravity duality and the mini-superspace approximation to the spectrum 
of certain hadronic states in the dual string theory. We then calculate the spectrum of eigenvalues and demonstrate that the level 
spacing is well described by the GOE; we also comment on an observation regarding the spectrum that might be a smoking gun 
of holographic models. We provide a series of appendices with technical details of various calculations. 



II. QUANTUM CHAOS IN THE HOLOGRAPHIC HADRONIC SPECTRUM 

The gauge/gravity correspondence has provided a very special window into the nature of strongly coupled field theories by 
mathematically identifying them with dual string theories fT0lfT2ll . The correspondence draws on a set of ideas dating back to 
't Hooft who argued that, in the appropriate limit, a field theory can be described by strings lfl3l . Since its inception, more than 
a decade ago, an important goal has been the extension of the correspondence to the case of theories with properties similar 
to those displayed by Quantum-Chromodynamics (QCD) in the strongly coupled regime; hoping to explain properties that are 
experimentally verified but defy a theoretical explanation, like confinement and chiral symmetry breaking. The gauge/gravity 
correspondence has produced various models of confinement and highlighted the universal behavior that follows from demanding 
that the dual Wilson loop satisfy the area law lfT4l . Remarkable progress in matching to the spectrum of light hadrons in QCD 
has been achieved using one of these models lfT5ll . 

An organizing principle of holographic dualities is that the energy scale in field theory is geometrized by an extra spatial 
dimension on the gravity side. Objects localized at large values of the r-coordinate are identified with states localized in the 
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ultraviolet of the field theory. Similarly, objects localized at small values of the r-coordinate are identified with states defined in 
the infrared of the field theory. Composite hadronic states that emerge from strong dynamics live mostly in the neighborhood of 
the smallest possible values of the r-coordinate, r m ^ n , which are naturally identified with the strong coupling scale Aqcd in the 
field theory. 




FIG. 1 : The holographic principle associates the direction r in the bulk with the energy scale in the field theory. 

In the gauge/gravity paradigm a classical trajectory in string theory serves as a description of a quantum state in the dual field 
theory. The oldest realization of this paradigm is given by Regge trajectories which precedes the AdS/CFT by more than thirty 
years. The old argument explaining Regge trajectories as spinning strings is now modified in the gauge/gravity correspondence 
(see Fig. ([]])). The classical trajectory whose conserved quantities (angular momentum J and energy E) satisfy J ~ E 2 is 
interpreted as describing the principal Regge trajectory, that is, a sector of the spectrum of hadrons with the same relation where 
now J is the spin and E — M is the mass of the hadron. Crucially, the string spins deep into the bulk geometry. The gravitational 
red shift between the small r region and the boundary (r — » oo) converts the string scale a' to the QCD scale. More precisely, 
the QCD string tension is given in terms of the geometry as Tqcd ~ #oo (jmm)/ (27rc/) highlighting the role of the holographic 
direction. 



A. Quantum spectrum from minisuperspace in string theory 

In standard quantum mechanics, starting with a classical picture there is a definite prescription for finding the energies of 
a system. We start with a classical Hamiltonian H, promote the generalized momenta and coordinates to operators and solve 
the Schrodinger equation. The eigenvalues are the allowed energies of the system. In our new, holographic picture, we are 
faced with new classical data - the trajectories of strings. How are these trajectories to be mapped to a quantum problem, and 
how do we extract information on the resulting spectrum? In string theory the Virasoro constraint provides the analog to the 
Schrodinger equation. It leads to the mass shell condition and is usually written as (H = (Lq — a)) \Sfr >= 0, where L is a 
Virasoro generator and a is a constant resulting from normal ordering. In this framework we think of > as a the precise 
analog of the wave function. 

We will follow a simplification known as the mini-superspace, an approximation where only modes corresponding to the 
center of mass are retained. The minisuperspace idea was originally introduced in the context of quantum cosmology fT6l ; in 
string theory, or more precisely, Liouville theory it was formulated in fFJj . More recently the mini-superspace has played a 
clarifying role in more complicated setups lfT8l . 
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We start with the Polyakov action in the conformal gauge: 

2 h I dT I da Gmn {~ xMxN + x ' Mx ' N ) > (d 



where M = 0, . . . , 9 and l/27ra' is the fundamental string tension. To implement the minisuperspace approximation, we 
consider an Ansatz where only one function depends on the world sheet coordinate a, let us call it X 9 and the rest depend 
on time, r: X M = (x n (r), X 9 (a)). After integrating equation |TJ) with respect to a we obtain essentially a particle-like 



Hamiltonian of the form 

H 



11 f 

= ^g^PmPn + ^V(x), g m n(x) / (1(tGmn(x), 



L 

\ '(./•) = J d(j(;^[.r)X n) x" ] . (2) 

o 



where p m are the canonical momenta conjugate to the coordinate x m and we have assumed that the metric Gmn is diagonal. 
The quantum spectrum is determined, in the framework of the minisuperspace formalism, by the following equation: 

- + V(x)V = 0, (3) 

where the Laplacian is computed in the metric gmn(%) given by equation This approach accesses only a special kind of 
string states (quantizing the center of mass motion). The minisuperspace also leaves out other important sectors of superstring 
theory like those coming from fermionic fields or fluxes in the geometry rendering the results somewhat limited. Nevertheless, 
in some specific situations with a large amount of symmetry the full answer of the string spectrum can be obtained from the 
minisuperspace approximation [19, 20l . 



B. Energy eigenvalues 

In what follows we will solve ^ for two supergravity backgrounds dual to confining field theories. 

The main idea of the holographic dictionary is to replace the study of a 4d field theory, in the appropriate limit, by that of a 
lOd gravity theory [10]. In this paper we are interested in lOd metrics of the form: 

ds 2 = A 2 (r)(-dt 2 + dR 2 + R 2 d(p 2 + dz 2 ) + B 2 (r)dr 2 + dsl (4) 

where r is the so called holographic direction, ds 2 represents a compact 5d space and can depend on r; its integrated volume form 
is denoted by u^(r). The space where the field theory lives (R 1 ' 3 ) has been parametrized in cylindrical coordinates (£, R, <p, z). 
The precise backgrounds we consider here are due to Maldacena-Nunez (MN) l2T1l and Witten (WQCD)[22| (see appendix |A| for 
details). 

Rather than a closed spinning string at the minimum of r corresponding to the holographic Regge trajectory (depicted in Fig. 
[TJ, we will consider a winding string which is largely localized in the same region, r m i n . Closed strings are generically mapped 
to glueball states; our spectrum will correspond, consequently, to the spectrum of glueballs. We hope that this spectrum forms 
the dominant part of the spectrum of hadrons but we are not able to directly apply our results to all hadrons. 

We focus on states for which (p(a) = a a, that is, where the string winds a times around the <p direction. A priori a is not a 
good quantum number in the backgrounds we consider as the string can unwind; we return to this issue later. In the formalism 
described in the previous section we set X 9 (a) — (f(cr). After some simplifications (see appendix [bJ we obtain a general 
problem of the form 



(-d 2 R - l(r)d r (m(r)d r ) + uo 2 R 2 A{rf) 9 = E 2 9, (5) 

where l(r) = l/(A(r)B(r)u^(r)), m(r) = A(r) 3 u^(r)/ B{r) and uo = a/(7ra / ). The functions A(r),B(r) anda^r) enter 
in the 10-d metric Eq. Q. To develop intuition into this spectral problem we comment on some of the most salient features at 
the heart of the results. The explicit form of the effective potential, V(p, R), in coordinates (p(r), R) where the kinetic terms 
are canonical, can be found in appendix [b] For both theories analyzed, the potential V(p, R) is impenetrable as p approaches 
its spatial boundary. Near p = the MN potential has a finite piece which is independent of R while the WQCD potential 
blows up as —1/p 2 . The point that we would like to emphasize graphically is that the potentials provide a rough realization 
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of a Bunimovich stadium l23l (see figure [2|. In the figures we exploit the symmetries of the problem around p = where 
in the full geometry a closing cycle leads to a smooth origin of the form R 3 for MN and R 2 for WQCD. With this symmetry 
configuration in mind we should characterize both potentials as Bunimovich-like stadia with a bump in the center, p = 0, a la 
Sinai. The similarity of our potentials with the typical potentials of quantum chaos (Bunimovich and Sinai) makes the appearance 
of quantum chaos more plausible. 




FIG. 2: Density plots of the effective potential, V(p, R), in the MN and WQCD backgrounds. 



III. THE SPECTRUM OF EXCITED HOLOGRAPHIC GLUEBALLS 

Quantum chaos studies the quantum properties of classically chaotic systems. In such systems the local statistics of the 
energy spectrum play a key role. An important feature is the level spacing distribution P(s), which is the distribution function 
of nearest-neighbor spacings E n +i — E n as we run over all levels plays a key role. A dramatic insight of quantum chaos is 
given by the universality conjecture for P(s). If the classical dynamics is chaotic, then P(s) coincides with the corresponding 
quantity for the eigenvalues of a suitable ensemble of random matrices P(s) = s exp(— s 2 /M 2 ) Q. 

We solve the mini-superspace spectrum ^ for two prominent supergravity backgrounds dual to confining theories: MN and 
WQCD (see appendix |A| for details of the backgrounds.). The eigenvalue problem itself is solved numerically using spectral 
decomposition [24]. We relegate most technical details to the appendices [B] but note a few key features of the problem. In the 
spectral method, the left hand side of Eq. ([5} is written as a 2-dim matrix folding the two spatial dimension R and p together 
using Kronecker outer product. With a slight abuse of notations we call that matrix the Hamiltonian. The energies are the 
square roots of the Hamiltonian eigenvalues. We exploit the parity symmetry in the R direction to sort the states by parity and 
conserve computational resources. Our basis for spectral decomposition in the indirection uses Hermite functions, since both 
d\ and R 2 appear as tridiagonal matrices. The complete Hamiltonian is non-zero only for elements that are within the N r + 1 
sub-diagonal (where N r is the number of spectra functions used in the r direction). The block-diagonal form of the Hamiltonian 
is very different from the dictum of Wigner-Dyson assumptions as summarized by RMT. We will track the implications of this 
difference. 

After using the reparametrization freedom in R and p in both problems we are left with a single parameter a/ (ira f ) which 
measures the interaction strength, we arbitrarily set a/ (na') = 1. This paremeter sets the energy scale at which the interaction 
term becomes important. Below that scale the Hamiltonian is, to a good approximation, a direct product and the eigenvalues will 
follow an uncorrelated Poisson distribution. 

In figure ^ we plot the absolute value square of a typical wave function for both backgrounds. The fast decay in p is the 
direct consequence of the form of the potential in that direction (compared to the quadratic behavior in the R direction). In 
both cases we can see that the value of the wave function near R = is small, this is not an accidental property of the Hermite 
functions. Rather, this indicates that the interaction term excites a large number of modes in the R directions. This observation is 
reassuring, the strings we consider carry a fixed winding number, denoted by a before, around the origin. This winding number 
is not conserved near the tip, R = 0, where the strings can unwind. Since the winding number is a good quantum number for 
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strings that are far away from the tip, we expect the mixing of different winding modes to be proportional to the density of the 
wave functions near the tip which is shown in fig. ([3} to be small. Therefore the wave function smallness near R = suggests 
that the mixing of different winding sectors will be a small effect. 



Energy = 14.4432 Energy = 16.3103 




FIG. 3: Wave functions. Smallness near R — guarantees that the mixing among states of different winding number is negligible. 

After calculating the spectrum for both cases we are now ready to study the eigenvalue spacing distributions. First we choose 
an arbitrary energy range (above the range of small interactions) spanning about 400 energy levels. We calculate the energy 
difference between eigenvalues in the chosen range and plot them on an histograms, and compare them to Wigner distributions 
P(s) ~ s exp( — (s/M) 2 ). The results are plotted in Fig. The root mean square (RMS) between Wigner's distribution and 
data is below 10 -3 when the distribution is normalized so the sum is one. This excellent matching to Wigner distribution proves 
our main claim that the spectrum of hadrons in the MN and WQCD theories shows a quantum chaotic eigenvalue distribution. 




FIG. 4: Eigenvalue spacing distribution P(s) for the MN background and for the WQCD background. 

It is worth noting that the lowest lying eigenvalues lead to a Poisson distribution and we drop them from consideration. This 
is not the case in QCD but the discrepancy can be attributed to the fact that holographic models require a parametrically large 

We repeat the above procedure for the entire energy spectrum we can accurately calculate. Varying energy ranges and bin sizes 
we determine the value of the parameter M in the Wigner distribution as a function of energy (the mean energy of the energy 
range used for the histogram). We keep only the cases where the matching to the Wigner distribution was good (measured by 
RMS < 10 -3 ). The results are displayed in Fig. §5§. The main feature shown is M decreasing with energy. This property 
differentiates between the models we study and the random matrices discussed in the Wigner-Dyson approach (where M does 
not depend on energy). We attribute this result to the clustered block-diagonal nature of the Hamiltonians, however we cannot 
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prove this connection. We conjecture that the dependence on energy M(E) is universal to all string theory models for confining 
theories beyond the two we study. 




FIG. 5: Dependence of M on the energy region for the MN background and for the WQCD background. 



IV. CONCLUSIONS 

As in previous studies, the holographic direction, r, has been crucial. Simple models of strings in flat space such as the original 
string approach of the 1960's would not lead to the results presented here. The condition that geometrically the supergravity 
backgrounds have to satisfy two properties which are sufficient for an area law for the dual Wilson loop has been our only input. 
Most of our results follow from the "end of the space wall", that is, the fact that goo {fmin) has a non- vanishing minimum at 
the smallest value of the holographic direction r. This behavior implies that for supergravity backgrounds like AdS$ x S 5 dual 
to conformal theories, the spectrum will be completely different. As shown in [25 ], the holographic area law condition alone 
implies chaotic behavior in the classical equations of motion and we have now explicitly verified that the quantum spectrum is 
chaotic in the standard sense of the GOE distribution. 

We have discussed two of the most studied models of holographic confining theories, the MN and WQCD models. We have 
found that the spectrum of nearest-neighbor eigenvalues has a Gaussian distribution consistent with the GOE. Interestingly, this 
is the distribution observed in nuclei, hadrons and realistic lattice studies of QCD and its supersymmetric versions. 

More speculatively, we have pointed out to a property that might be specific to holographic models of confinement: the 
dependence of M on the energy range E, in the statistics of eigenvalue levels. 

We have thus, using the principles of holography, provided a concrete Hamiltonian for hadronic excitations with the properties 
conjectured by Wigner more than half a century ago. In a a holographic sense we have provided a concrete realization of Wigner's 
ideas. This is thus an explanation of quantum chaos or RMT in a sector of the spectrum of the hadrons in holographic models. 
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Appendix A: Supergravity backgrounds 

1. The Maldacena-Nunez background 

The MN background whose IR regime is associated with J\f = 1 S YM theory is that of a large number of D5 branes wrapping 
an S 2 ['21 ] (see also (261). To be more precise: (i) the dual field theory to this SUGRA background is the J\f = 1 SYM 
contaminated with Kaluza- Klein (KK) modes which cannot be de-coupled from the IR dynamics, (ii) the IR regime is described 
by the SUGRA in the vicinity of the origin where the S 2 shrinks to zero size. The full MN SUGRA background includes 
the metric, the dilaton and the RR three-form. It can also be interpreted as an uplifting to ten dimensions a solution of seven 
dimensional gauged supergravity [26 1. The metric and dilaton of the background are 

ds 2 = [r]^dx^dx v + a'g s N (dr 2 + ds 2 )] , 

ds 2 =e 2 ^\eUe 2 )^\(e 2 ^eUe 2 )^ 
sinh 2r 



e 20 = e -20 o 



2e#( r ) ' 

r 2 1 



where /i, v = 0, 1, 2, 3 and 



e 2s(r) =rcot h2r , , (Al) 

sinh 2 2r 4' 



ei =d6i, 62 = sin#i<i</>i, 

es = cos ip d$2 + sin ip sin 62 d(f)2 — a{r)d6\ , 

e4 = — sin ij; d02 + cos -0 sin 62 d(j)2 — a(r) sin 0id(j)i , 

=d-0 + cos 62 d<j)2 — cos a(r) = — r 2 , (A2) 

sinh r 



where fi = 0, 1, 2, 3, we set the integration constant e^ D o = y/g s N The 3-form is 

1 



ff 1 ^ = g s N 



^(w 1 — A 1 ) A (w 2 — A 2 ) A (w 3 — A 3 )- 



a 

A = ^ [cr 1 a(r)dfli + a 2 a(r) sin 0id0i + a 3 cos (A3) 

and the one-forms w a are given by: 

w 1 + uu 2 = e~ irtp (d6 2 + zsin0 2 #2) , = # + cos0 2 #2 (A4) 

Note that we use notation where have dimension of length whereas p and the angles #i, </>i, #2? 02? ^ are dimensionless 

and hence the appearance of the a 7 in front of the transverse part of the metric. 

2. The Witten QCD background 

In this section we present the Witten QCD background; we include some notes on the holographic relations with the field 
theory (22). The ten-dimensional string frame metric and dilaton of the Witten model are given by 



ds 2 



3 

f(r)= 1- j|, L = (7rNg s )ia' i , 

e* = g s (y • (A5) 
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The geometry consists of a flat, warped 4d part, a radial direction r, a circle parameterized by 6 with radius vanishing at r = rrj, 
and a four-sphere whose volume is instead everywhere non-zero. It is non-singular at r = 7*0. Notice that in the r —> 00 limit 
the dilaton diverges: this implies that in this limit the completion of the present IIA model has to be found in M-theory. The 
background is completed by a constant four-form field strength 

F 4 = 3L 3 uj 4 , (A6) 

where CJ4 is the volume form of the transverse S 4 . 

An important gauge theory parameter is the KK mass scale I /Re, which is given by 

= 2 m ° ' where m o = ^ • ( A7 ) 

As can be read from the metric, mo is also the typical glueball mass scale and, its square is proportional to the ratio between 
the confining string tension Tqcd and the UV 't Hooft coupling A. As usual, the supergravity approximation is reliable in 
the regime opposite to that in which the KK degrees of freedom can decouple from the low energy dynamics. The condition 
Tqcd <^ ^0 implies in fact Ac 1, which is beyond the supergravity regime of validity. 

Appendix B: Mini-superspace quantization of winding strings 

In this section we elaborate on the details of the mini-superspace quantization of the WQCD and MN theories. First we rewrite 
the flat subspace in ( |A1| ) and ( |A5| ) using cylindrical coordinates 

r]^ v dx^dx v = -dt 2 + dR 2 + R 2 dtp 2 + dz 2 (Bl) 

Singling out (p as the direction the string will wind around, the mini-superspace quantization, Eq. ([3]), reduces to 

• MN-background 



2 

E 2 - -e 2 * R 2 ) ¥ = (B2) 



• WQCD-background 

+ j^rr^ ( rl9/4 /( r ) + d R*+ 

E - - — R 2 ) * = 0. (B3) 



(ira>) 2 L 3 



We can 'symmetrize' the equation using the following redefinition of the wave function 



*(R,r) = -±-MR,r) (B4) 
r \r) 

with F{r) — ei^ g and F(r) = L^r^ fi (r) for the MN and WQCD theories respectively. This trick simplifies the Kinetic 
term in r however now the entire equation has an overall factor of F(r). This factor can be easily removed, however we now 
need to add a boundary condition of ip, such that it vanishes at zeros of F(r). In both MN and WQCD cases the wave function 
equation can be written in the following form : 

°~ k (r)-0-(k(r)4-)+p(r)+ 



F(r) 



d 2 R dr V dr 



+ (^-) 2 R 2 q (r)-E 2 



i/;(R,r) = 0. (B5) 



With 
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• MN 

k(r) =1 

p(r) = -^"{r)+g"{r) + Q^'(r) + </(r)) 

q(r) = e 20(r) (B6) 



• WQCD 

k(r) 



1 



L 3 / 2 



'o 



P(r)= L3 (-6^ (fc(r))2 + ¥ fc(r)fc ' (r)+ 
+ J(fc'W) 2 + ^(r)fc"(r)j 

3 

= ^ (B7) 



The transformation to the canonical coordinate used in Fig. ^ is 

dv 

= dp, (B8) 

Moving to the canonical coordinates is intuitive but not necessarily advantageous from a numerical point of view. In the canonical 
coordinates the large r behavior is very steep and standard orthogonal polynomials will not give good convergence. Instead we 
use different set of coordinate that allow us to use Hermite polynomials: 

• MN 

dr dp 

W) ~ T^+F 

• WQCD 

dr dp 



k(r) / p 4 +6p 2 r + 12^ 

8L 3 



Appendix C: Finding the energy spectrum using Spectral methods 



In this appendix we describe the solution eigenvalue problems ( |B5| ) using spectral methods. For a detailed description of 
spectral methods see l24l . First using reparametrization we set (j) = g s N = 1 and L = ro = 1 in MN and WQCD, reducing 
both cases to a single parameter w = a/ (aV). We expand the wave-functions using the wavefunctions of harmonic oscillator 

R\ , I p 



iy y / 
T )= n m— n ^ 



n=0 m=0 

1 Z? 2 , 



7 / V 7 



ip n (x) = e i R h n (x) , 

^L m {x) = e"^ 2 hz m (x), (CI) 

where h m (x) are the normalized Hermite polynomial and hz m (x) are a modified version of Hermite such that hz m (0) = 0. 
Both sets of orthogonal polynomial can be calculated from a recursion relation: 



h n +i(x) =\l — ~~r xh n (x) + a I — ^— h n _i(x), (C2) 



lE 2 7i+2(^) =\/ 2n ^ +3 ^hz2n+l(^) + y 2 n + 3 hz 2n (^), 
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with the initial elements 



hzi(x) = hi(x), 

Solving ( |B5| ) reduces to finding eigenvalues and eigenvector of the matrix 



h\{x) = V2xho(x), 
2 



hz 2 (x) 



xhzi(x) 



(nm)(n'm' 



_ f 



dRdp , ( R\ 

R2 I Rip \ L R, 



-p(p) Wirg(p) 



(C3) 



In the paper we use results for w = 1 after checking that a order of magnitude change in w give qualitatively similar result. We 
calculate the energy spectrum for several values of grid size Nr = N p = 20, 40, 60, ... , 160. This allow us to evaluate the 
error of the calculation by comparing result with different grid sizes. The choice of Ir and l p is done by optimizing for minimal 
errors. When the scales are properly chosen, the errors are proportional to Nr * N p which matched the expectation. At large 
energy the wavefunction gets nearer to the large r region, where the expansion is bound to fail. The reason behind that is that the 
equation is not separable even at large r so it is impossible to get a good expansion based on a rectangular grid of 1 -dimensional 
functions. The effect of this is that the number of reliable eigenvalues does not grow as Nr * N p but rather saturates (in our case) 
near 1500 for MN and 1000 for WQCD. Luckily the number of eigenvalues we reliably obtain is large enough for the statistical 
measures described in the main text. 
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